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Abstract 
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nonzero chemical potential. These effects are important for cosmological studies of the early Universe. 
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I. INTRODUCTION 



The quark-gluon plasma (QGP) is a state of matter which forms at very high temperatures 
or densities. It is believed that up to microseconds after the big bang the QGP was a dominant 
component of the Universe. This state of matter is recreated in heavy-ion collision experiments 
[such as are done at the Relativistic Heavy Ion Collider (RHIC)] which study its formation and 
transition to ordinary matter. The equation of state (EOS) of the QGP is essential to our under- 
standing of its hydrodynamic expansion and consequently of the particle spectra produced in these 
experiments. We have studied the EOS at zero and nonzero chemical potential previously yj, |2fl. 
Here we extend our work in two directions. (1) We present results for the EOS at nonzero chem- 
ical potential at finer lattice spacings than our previous work. Here the temporal lattice extent is 
N t = 6, where previously it was N t = 4. Preliminary results for the N t = 6 case were reported in 
Ref. [3]. It is important to compare the two cases and determine the size of the discretization error 
as a step towards taking the continuum extrapolation. (2) We include the effects of the charm and 
bottom quarks. A preliminary progress report on the charm quark effects was given in Ref. |4J]. 
We use the heavy-quark-quenched approximation. That is, the charm and bottom quarks appear 
as valence quarks, but not as dynamical sea quarks. Thus we ignore all charm and bottom quark 
loops contributing to the operators we determine in order to obtain the EOS. This approximation 
introduces an error in our calculation. However, considering that the charm and bottom quarks are 
much heavier than our sea u, d, and s quarks, it seems plausible that adding sea charm and bottom 
quarks would have a small effect for temperatures much less than their masses. Still, until we have 
a dynamical c- and ^-quark calculation to compare against, this statement remains a conjecture. 
The equation of state with the charm and bottom quarks added is most applicable to the study of 
the early Universe, since the time scale relevant to the heavy-ion collisions at RHIC is probably 
too short for the charm and bottom quarks to thermalize and have a visible effect on the particle 
data. 

As in our previous N t = 4 determination of the EOS at nonzero chemical potential, we employ 
the Taylor expansion method. For a detailed description of the method, see Refs. [5, (]]. The 
expansion is carried up to sixth order in the expansion parameters /u q /T, where /u q is the chemical 
potential for a certain quark flavor q and T is the temperature. 

The gauge ensembles we used in this work are the same as in Ref. Q2Q . They are generated 
using the asqtad improved staggered action [|7D and have two degenerate light quarks and a strange 
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quark in the sea. The ensembles lie approximately on a trajectory of constant physics, where the 
strange quark mass m s is tuned to be close to its physical value, and the light quark mass m/ is 
one-tenth of m s . Because in this paper we also consider charm and bottom quarks, we do not refer 
to the strange quark as the "heavy quark" as in Refs. Jl, 2]. 

In Sec. II, we present our results for the 2+1 flavor EOS with nonzero chemical potential at 
Nt = 6, and compare it with our previous one at N t = 4. We also show other thermodynamic 
quantities, such as the quark number susceptibilities and light-quark density. Section III gives 
our findings for the isentropic EOS for 2+1 flavors. In Sec. IV, we calculate the effects of the 
charm quark on the EOS at zero and nonzero chemical potential, using the heavy-quark-quenched 
approximation to represent it. Section V does the same for the bottom quark. In Sec. VI, we give 
our conclusions. The Appendixes contains some helpful formulas for the application of the Taylor 
expansion method for the EOS calculation in the 2+1+1 quark flavor case. 

II. THE EOS AT NONZERO CHEMICAL POTENTIAL AT N, = 6 FOR 2+1 FLAVORS 

The Taylor expansion method allows us to represent the pressure p and the interaction measure 
/ in the case where both the light and the strange quark chemical potentials are nonzero, as the 
following infinite sums: 

P _ lnz _ ^ (r\(^L\ n (VlC m 

T 4 - T 3y ~ 1* CnmK 1 ) y T J ^ J , 

In the above, are the chemical potentials for the light and strange quarks in physical units, T 
is the temperature, N s is the spatial lattice extent, Z is the partition function and a is the lattice 
spacing. Because of CP symmetry, the expansion coefficients c nm and b nm are nonzero only if 
n + m is an even integer. The explicit forms of these coefficients are given in Ref. 020, and since they 
are somewhat involved, we do not repeat them here. We calculate the coefficients stochastically 
with random Gaussian sources. Inside the transition region, we used 800 random sources per 
lattice, and outside, 400. With these numbers, the stochastic error in the unmixed second order 
coefficients, i.e., the diagonal quark number susceptibilities at zero chemical potential, is about 
20% of the full statistical error. These coefficients are the ones with largest contribution to the 
thermodynamic quantities for each type of quark. For the fourth order unmixed coefficients, this 
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FIG. 1 : Unmixed coefficients c„q and cq„ in the Taylor expansion of the pressure as a function of tempera- 
ture. The new results for N t = 6 are shown in filled (red) circles; empty (black) circles are used for N t =4 
(from Ref. I2D). Arrows indicate the Stefan-Boltzmann limit for each of the coefficients. 

error is about 50% of the final error. For the rest of the coefficients (mixed coefficients and all 
coefficients of sixth order), the contribution of the stochastic error is dominant. A further increase 
of the number of sources as a way to decrease the stochastic noise seems impractical at this point. 
We need either significantly more computer power or a substantial improvement of the noisy 
estimators in order to reduce the resulting stochastic error. 

Figures [T] and |2] show some of the coefficients in the pressure expansion and compare our new 
results at N t = 6 (red filled circles) with the previous ones [2] at N t = 4 (black empty circles). We 
can see that the errors for the Nt = 6 case are smaller than the ones at the shorter temporal extent, 
due to both the increased volume and increased number of random sources. (We previously used 
100-200 sources.) There is also a shift in the central values between the two cases which indicates 
that the discretization effects at Nt = 4 are significant. The approach of the coefficients to the 
(massless) Stefan-Boltzmann continuum limit with increasing T in the case of N t = 6 is slower, 
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FIG. 2: Mixed coefficients c nm in the Taylor expansion of the pressure as a function of temperature. The 
new results for N t = 6 are shown in filled (red) circles; empty (black) circles are used for N t = 4 (from 
Ref. [2]). Arrows indicate the Stefan-Boltzmann limit for each of the coefficients. 



and the structure at low temperature is made somewhat clearer due to the smaller errors on the data. 
Similarly, Figs. [3] and |4] compare the unmixed and mixed coefficients b nm involved in the Taylor 
expansion of the interaction measure at the two different temporal extents. These coefficients 
are calculated independently from the ones in the pressure expansion (although the two sets of 
coefficients are technically related by integration). The operators involved in the determination 
of the interaction measure expansion coefficients are intrinsically noisier, which is reflected in the 
larger errors on the data. Still, the shift in central values and errors in the b nm coefficients in going 
from N t = 4 to 6 is qualitatively similar to that of the pressure coefficients c nm . 

Having obtained the coefficients in the Taylor expansions in Eqs. (OQ) and ©, we can now turn 
to calculating the EOS. Because of the nonzero c n \ (T) terms, a nonzero strange quark density n s is 
induced even with /j s = 0. (We use /jf to denote the chemical potential in lattice units for flavor /). 
To study the n s = plasma, we must, therefore, tune ju s as a function of /ji and T. Figures \5\ (both 
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FIG. 3: Unmixed coefficients b n o and bo„ in the Taylor expansion of the interaction measure as a function 
of temperature. The new results for N t = 6 are shown in filled (red) circles; empty (black) circles are used 
for N t = 4 (from Ref. Arrows indicate the Stefan-Boltzmann limit for each of the coefficients. 

panels) and [6] (left panel) show the changes in the interaction measure (A/), pressure (Ap) and 
energy density (Ae — AI + 3Ap) at fii/T = 0.1, 0.2, 0.4, and 0.6, and with fi s /T tuned along the 
trajectory so that n s rs 0. We see statistically significant discretization effects when we compare 
the Nt = 4 and 6 cases for Ap, with the latter data lying lower than the former. The AI results have 
larger errors, and discerning differences in the data from the two temporal extents is more difficult. 
The N t = 6 data is slightly but consistently lower than the one from the N t = A calculation; however, 
this is not a statistically significant observation. The change in the energy density Ae inherits the 
large errors from AI, and the same conclusions apply for it. The discretization effects for the light- 
quark density (n uc {), the light-light and strange-strange quark number susceptibilities (% J([ , and % ss , 
respectively), are examined in Figs. [6] (right panel) and [7] (both panels). The effect of increasing 
the temporal extent from N t = 4 to 6 is to lower these quantities by 4%-10%. For the range of 
values of jji/T that we examine, we do not find any evidence for peaks that could presage critical 
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FIG. 4: Mixed coefficients b mn in the Taylor expansion of the interaction measure as a function of tempera- 
ture. New results for N t = 6 are shown in filled (red) circles, empty (black) circles are used for N t =4 (from 
Ref. 11). Arrows indicate the Stefan-Boltzmann limit for each of the coefficients. 



behavior in % uu . The light-strange quark number susceptibility (% t(5 ), shown in Fig.[8](left panel), is 
too noisy for a reliable conclusion about its discretization effects; there is only a hint at a possible 
move toward lower absolute values at the larger N t . 



III. THE ISENTROPIC EQUATION OF STATE 

The form of the EOS most applicable to the experimental conditions of the heavy-ion collisions 
is the isentropic one. There, after thermalization, the system expands and cools with constant 
entropy. We determine the isentropic EOS by performing our calculations at fixed ratio of entropy 
to baryon number (s/hb). This is achieved by finding the trajectories in the Qui, /u s , T) space which 
satisfy (within errors) both s/ng = C and n s = 0, where C is a constant whose value depends on 
the particular experiment we are interested in. For AGS, SPS, and RHIC, we have s/ng = 30, 
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FIG. 5: (Left panel) The change in the interaction measure due to the nonzero chemical potentials vs. 
temperature. At a given N t , the larger fii/T, the higher the data appears on the plot. (Right panel) Similarly, 
the change in the pressure. 
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FIG. 6: (Left panel) The change in the energy density due to the nonzero chemical potentials vs. tempera- 
ture. (Right panel) The light-quark density vs. temperature for several values of jiii/T (in different colors) 
and fj s tuned such that n s « 0. 

45, and 300, respectively. Figures [8] (right panel) and [9] (both panels) show our results for the 
interaction measure, pressure, and energy density for the different s/hb values appropriate for 
these experiments. We compare the new results at Nt = 6 (filled symbols) with the ones already 
published in Ref. [2] at N t = 4 (empty symbols). For these three quantities, the comparison shows 
negligible effects due to the increase of the temporal extent N t . The reason for this is that by 
far the largest contribution to these quantities — the zero-chemical potential (zeroth order) term in 
their respective Taylor expansions — does not show large discretization effects [ l|. On the other 
hand, quantities which do not have a zeroth order term may show larger differences between the 
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FIG. 7: (Left panel) The light-light quark number susceptibility vs. temperature for several values of 
(in different colors) and /u s tuned such that n s ns 0. At a given N t , the larger fii/T, the higher the data 
appears on the plot. However, at low fii/T these differences are very small. (Right panel) Same for the 
strange-strange quark number susceptibility. Here the data dependence on jui/T is very weak. 
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FIG. 8: (Left panel) The light-strange quark number susceptibility vs. temperature for several values of 
jui/T (in different colors) and /u s tuned such that n s « 0. As in Fig. [7] (right panel), the data dependence on 
fii/T is very weak. (Right panel) The isentropic interaction measure vs. temperature for selected values of 
s/n B . 

N t = 6 and 4 cases. Indeed, small discretization effects are evident in the isentropic light-light 
and strange- strange quark number susceptibilities in Figs. [10] (right panel) and [TT] (left panel), 
respectively. However, the isentropic light-quark density, shown in Fig. [101 (left panel), has only 
marginal discretization effects, despite the fact that it does not have a zeroth order contribution. 
The large errors on the strange-light susceptibility, shown in Fig. [TT] (right panel), precludes us 
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FIG. 9: (Left panel) The isentropic pressure vs. temperature for selected values of s/ng. (Right panel) The 
same for the isentropic energy density. 
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FIG. 10: (Left panel) The isentropic light-quark density vs. temperature for selected values of s/ub- (Right 
panel) The same for the isentropic light-light quark number susceptibility. 



from drawing a conclusion about its discretization effect. To conclude with our final observation, 
for both values of N t , we find rather smooth behavior for the isentropic variables indicating that 
experiments are far from any critical point in the p — T plane. 



IV. THE EFFECTS OF THE CHARM QUARK ON THE EOS 

In this section, we study the effects of the charm quark on the EOS at zero and nonzero chemical 
potential. Our preliminary results were reported in Ref. [J]. First, let us discuss the relevance of 
the charm quark contribution. The experiments at RHIC create a "fireball" which thermalizes 
within x w 10~ 24 s [[8|]. The u, d, and s quarks participate in the thermal ensemble describing the 
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FIG. 1 1 : (Left panel) The isentropic strange-strange quark number susceptibility vs. temperature for se- 
lected values of s/tis- (Right panel) The same for the isentropic light-strange quark number susceptibility. 

state of the thermalized fireball. Under the experimental conditions, the c quark probably is not 
thermalized, and thus the 2+1 flavor EOS is considered sufficient for the hydrodynamics models 
applied to the current experimental data. The question of equilibration of charm, however, is not 
completely settled as argued, for example, in Ref. Q9[|. Furthermore, the situation may change for 
the future LHC experiments. A quark-gluon plasma also existed microseconds after the big bang. 
Under these primordial conditions and longer time scales, the c quark probably participated in the 
thermal ensemble as well, which implies that for the study of the early Universe, the EOS with 
2+1+1 flavors would be important dlOfl. For example, the scale factor of the early Universe is 
affected by the number of quark flavors in the EOS used for its determination [11]. Previously, the 
question of the charm quark contribution to the EOS at zero chemical potential has been studied 
on the lattice in Ref. [12] at N t = 4, 6, and 8 using the p4 fermion formulation. That study treated 
the charm quark as a valence staggered quark. We do the same, but in the asqtad formulation at 
Nt = 6. We tuned the charm quark using a different strategy than in Ref. [|12|l . where the charm 
quark mass was determined using the T] c . or J /m rest mass on all available ensembles. In our study, 
the charm quark mass was tuned to match the rest mass of the D s at (3 = 7.08 {a ~ 0.086 fm) where 
the discretization effects are smallest on our trajectory. We chose the D s for our tuning purposes 
because the discretization effects are smaller for heavy-light mesons than for the heavy-heavy 
ones |13|1. We found m c /m s = 10 at our tuning point with a 4% uncertainty. We have kept this 
ratio constant for lower temperatures. It is probably incorrect at the lowest available temperatures, 
but due to large discretization effects, the tuning is inherently problematic there. Still, we do not 
expect this to matter much, thanks to the large mass of the c quark and its very small contribution 
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in that region. 

Following our method in Ref. [1] for determining the EOS at zero chemical potential, the 2+1+1 
flavor interaction measure was obtained by adding to our previous results for the 2+1 flavor case 
the charm contribution 



d(m c a) , ._ . dun . I dM 



(3) 



dlna dlna \ dun 

where the observables in the above are calculated in the heavy-quark-quenched approximation and 
A stands for the difference between the zero and nonzero temperature value of an observable. The 
mass beta function is approximated as 

d(m c a) _ ^ d(m s a) 
dlna dlna 

since we kept the ratio m c /m s = 10 constant along the trajectory. We determined the strange quark 
mass beta function and the function dun I dlna previously yj]. To find the charm contribution to 
the pressure and energy density, we integrated Eq. §3$ along the physics trajectory, as in Ref. yj] 
for the 2+1 flavor case. 

Again, the nonzero chemical potential calculation was done using the Taylor expansion method, 
taken to sixth order. For 2+1+1 quark flavors, the Taylor expansion of the pressure is modified to 
the following form: 

oo /-\n/-\m/-\k 

where jli jSjC are the chemical potentials in physical units for the light (u,d), strange (s) and charm 
(c) quarks. Because of CP symmetry the terms in the above are nonzero only if n + m + k is even. 
The interaction measure has the same form with only c nm t — > b nm k. Some details of the explicit 
calculations for the pressure and interaction measure coefficients can be found in the Appendixes. 
We used 800 random sources per lattice in the transition region and 400 outside it to calculate the 
new observables in the expansions of the pressure and interaction measure. For the calculation at 
nonzero chemical potential, the valence c quark had a low cost in terms of computer time, but it 
required a sizable software development. For 2 + 1 flavors we had 95 observables to code and for 
2 + 1 + 1 flavors there were 399. 

Turning to our results, let us first examine the effects of the charm quark on the EOS at zero 
chemical potential. Figure [T2l (left panel) shows our results for the EOS with 2+1+1 flavors and 
compares it with previous results for 2+1 flavors [1]. The charm quark contribution grows with 
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FIG. 12: (Left panel) Interaction measure (/), pressure (p) and energy density (s ) divided by the temperature 
to the fourth power (T 4 ) for the cases of 2+1 (red) and 2+1+1 (black) flavors. The arrows indicate the energy 
density Stefan-Boltzmann limit for both cases. (Right panel) The pressure for 1 quark flavor in the free 
theory vs. the ratio of the quark mass (m) and temperature (T) for different staggered quark formulations. 
The rise of the pressure at large m/T in the HISQ case shows that higher order corrections to the Naik term 
are needed in this region. Currently we have corrections up to 0(m & ) only. 

temperature, as expected, and at the highest available T it contributes about 20% to the energy 
density. We conclude that in the cases where the charm quark is thermalized, its contribution to 
the EOS at temperatures higher than about 200 MeV, cannot be ignored. Our result at N = 6 is 

n 

qualitatively similar to the previous work H12H . but quantitatively our charm quark contributions to 
the energy density and pressure are about 25%-30% lower by comparison at temperatures around 
400 MeV. A possible explanation for this is the larger discretization effects for the heavy-quark 
pressure for the asqtad action than for the p4 action. Figure [12] (right panel) shows the free quark 
pressure as a function of the ratio of the (heavy) quark mass and the temperature for different 
staggered lattice fermion formulations. The asqtad action at Nt = 6 shows a negative value for the 
pressure for a range of heavy-quark masses while the p4 action is close to the continuum limit. 
Our results for the charm contribution to the EOS do not show the outright unphysical behavior 
occurring in the free quark case, but it is possible that the heavy-quark discretization effects depress 
the lattice values. 

Now let us turn to the results at nonzero chemical potential. Figures [13] and [14] present some 
of the pressure and interaction measure expansion coefficients which are directly related to the 
charm quark contribution at nonzero chemical potential. The first row in both figures shows the 
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unmixed coefficients and the second row — three of the mixed coefficients. The mixed coefficients 
are quite small and are much noisier than the unmixed ones, which was expected. As a whole, 
the new unmixed coefficients cqq„ and Z?oon m me pressure and interaction measure expansions are 
small compared with the c„oo> c om, b n oo, and bo n Q coefficients. For numerical comparisons see 
Sec. II, where the latter four sets are defined without the last zero in the subscripts. These new 
coefficients remain well below the continuum (massless) Stefan-Boltzmann values at the highest 
temperature available here. This is not surprising, since over our temperature range T < 2T C , the 
charm quark mass is much larger than the temperature. The first panel of Fig. [13] shows that 
C002 becomes slightly negative for temperatures up to about 220 MeV. This behavior is obviously 
unphysical, since this coefficient is directly proportional to the necessarily positive charm quark 
number susceptibility at zero chemical potential % C c(w,j,c = 0) ~ ( n c)> where n c is the charm 
quark number density. It follows that C002 should be a non-negative number at all temperatures. 
We tracked this unphysical behavior to the interplay between the heavy-quark mass and the Naik 
term in the asqtad action. In the tuning of the latter, corrections proportional to were not 
included. It is easiest to understand this if we examine the quark number susceptibility for free 
asqtad (Naik) fermions at large quark masses shown in Fig. [TBI (left panel). In the continuum limit, 
this susceptibility should approach zero from above with increasing heavy-quark mass. We find 
that at Nt = 6 and 8 there is a pronounced "dip" into negative values for a certain range of large 
quark masses. This effect is much smaller at N t = 12. Since this particular discretization effect 
does not occur for standard staggered fermions at N t = 6, we conclude that certain thermodynamic 
quantities, such as susceptibilities, are sensitive to the "length" of the Naik term and require large 
Nt's in order to overcome their unphysical behavior. From Fig. [15] (left panel), the p4 action 
seems to be much closer to the continuum limit at N t = 6 and very probably will not show this 
particular discretization effect in the dynamical case. The HISQ action [14] improves the heavy- 
quark dispersion relation by tuning the coefficient of the Naik term. (The same tuning could have 
been done with the asqtad action.) Tuning suppresses this unphysical behavior for N t > 6 for the 
range of m/T up to 0(8). Still, in our unquenched 2+1 flavor case, the negative dip in the C002 
coefficient is quite small, so that its effect, for example, on the isentropic EOS is negligible over the 
parameter range relevant to heavy-ion collisions. Of course, other mixed and unmixed coefficients 
might be affected by the limited temporal extent N t = 6 as well, but since they are even smaller 
than C002 we can also ignore their unphysical contribution at low temperatures and small chemical 
potentials. 
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FIG. 13: Some of the new Taylor expansion coefficients for the pressure at nonzero chemical potential when 
the charm quark is added to the partition function. 

From the point of view of the isentropic EOS, our results show that the effect of the charm 
quark cannot be simply ignored. We have determined the approximate isentropic trajectories in 
the (jU/, fi s , fi c , T) space, by numerically solving the system 

S /2 c fir 

— (Ml,H s ,iu c ,T)=C, —r(fii,/u s ,p c ,T) = 0, —r(jui,/u s ,iu c ,T) = 0, (6) 

tiB t 3 r j 

with C = 30, 45, and 300. Figures [TBI (right panel) and [16] (both panels) present the 2+1+1 flavor 
isentropic interaction measure, pressure, and energy density, respectively, and compare them with 
the 2+1 flavor case. We see that the charm quark contribution is non-negligible, although it is due 
mainly to the contribution of the zeroth order coefficients in the Taylor expansions (i.e., the EOS 
calculated at zero chemical potential). We also note that for the range of temperatures between 
about 220 and 280 MeV, the errors on the isentropic interaction measure become large. In this 
region of the isentropic trajectory, /u c is big enough to make contributions from the quite noisy 
mixed coefficients visible. The isentropic energy density, of course, inherits this feature, being a 
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FIG. 14: Some of the new Taylor expansion coefficients for the interaction measure at nonzero chemical 
potential when the charm quark is added to the partition function. 

linear combination of the interaction measure and the pressure. 



V. EFFECTS OF THE BOTTOM QUARK ON THE EOS 

In the previous section, we presented evidence that the charm quark contributions to the EOS 
are non-negligible. At still higher temperatures the b- and eventually ?-quark contributions should 
be similarly non-negligible. In this section, we examine the effects of the bottom quark on the 
EOS in the range of temperatures up to about 400 MeV. Since the bottom quark is considerably 
heavier than the charm quark, we expect its contribution to the EOS to be smaller. To estimate 
it, we simply repeated the charm quark calculation but with a heavier mass corresponding to the 
bottom quark. The quenching error, even if relevant, will be smaller than the corresponding one 
for the charm quark. On the ensemble that we used for the charm quark, we tune the bottom quark 



16 



15 




N = 6, asqtad 




black: 2+1 + 1 
red: 2+1 



0.8 



■ — ■ N = 8, asqtad 
x-x N = 12, asqtad 
«~»N =6, p4 



0.6 



0.4 



N = 6, staggered 
►— ► N = 6, hisq 
— continuum 



0.2 








150 



200 



250 



300 



350 



400 



5 10 15 

m/T 



T [MeV] 



FIG. 15: (Left panel) The quark number susceptibility for 1 quark flavor in the free theory vs. the ratio of the 
quark mass (m) and temperature (T) for different staggered quark formulations. The reason for the rise of 
the susceptibility in the HISQ case for large m/T is the same as explained in the caption of the right panel of 
Fig- El (Right panel) The isentropic interaction measure at selected s/ns values for 2+1 and 2+1+1 flavors 
(red and black respectively). For a data set with the same color {i.e., produced with the same number of 
quark flavors), the highest lying results are for s/hb = 30, in the middle is the s/tib = 45 case and the case 
of s/hb = 300 has the lowest lying values. 

mass to match the B s rest mass to its experimental value. We found that within 5% mt>/m s = 
38. We kept that ratio constant along the physics trajectory. The problems of the tuning of the 
bottom quark are potentially worse than in the case of the charm quark, but we do not expect 
them to skew significantly our final result for the EOS, since the bottom quark contribution itself 
is expected to be small. Figure [FT] (left panel) shows the pressure and energy density at zero 
chemical potential with (2+1 + 1-1-1) and without (2+1 + 1) the bottom quark. We can conclude 
that the bottom quark contribution to the EOS at zero chemical potential is small (less than a 
standard deviation) in the transition region. It grows to about a standard deviation at temperatures 
close to 400 MeV. However, the range of temperatures we examine here is somewhat limited and 
probably by T ~ 600 MeV the bottom quark effects would grow to be statistically significant 
for comparable statistics at that temperature. We also have to bear in mind that the heavy-quark 
discretization effects may play a significant role here and keep the bottom quark contribution lower 
than what it would be in the continuum limit. 

As for the EOS at small nonzero chemical potential, our results for the coefficients of the Taylor 
series beyond the zeroth order term discussed above show that the bottom quark contribution can 
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FIG. 16: (Left panel) The isentropic pressure at selected s/ns values for 2+1 and 2+1+1 flavors (red and 
black respectively). The data ordering is as in the right pannel of Fig. Q3] (Right panel) The same for the 
isentropic energy density. 
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FIG. 17: (Left panel) The pressure and energy density vs. temperature for 2+1+1 and 2+1+1+1 flavors (red 
and black, respectively). The arrows indicate the energy density Stefan-Boltzmann values for both cases. 
(Right panel) The pressure Taylor expansion coefficients C0020 and C0002 vs. temperature. 

be safely ignored at the present level of statistics. The discretization effect which we found for the 
charm quark in the previous section is much worse for the bottom quark. Figure [FT] (right panel) 
compares the C0020 coefficient (referred to as C002 in the previous section) and the coefficient C0002 
in the Taylor expansion for the pressure when all chemical potentials /Ji ;S ,c,b 0- The C0002 is 
persistently negative at all available temperatures. At the present level of statistics and small 
nonzero chemical potential, this effect is unimportant. But if higher precision is desired, one 
should tune the coefficient of the Naik term and increase N t . 
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VI. CONCLUSIONS 



We extended our thermodynamics study of the quark-gluon plasma with chemical potential to 
finer lattices with temporal extent N t = 6. Comparing our results with previous results at N t = 4 
gives an indication of the importance of cutoff effects. As before, we used the Taylor expansion 
method to sixth order for the case of 2+1 quark flavors. We found small but significant changes 
in the coefficients of the Taylor expansions of the pressure and interaction measure in going from 
N t = 4 to 6. This leads to small differences in the resulting interaction measure, pressure, and 
energy density between the two cases, when matching to the experimental condition of zero strange 
quark density and keeping fli/T constant. Under these conditions, small discretization effects are 
also visible in the light-light, light-strange and strange- strange quark number susceptibilities and 
the light-quark density. On the other hand, the isentropic EOS shows very little difference between 
N t = 4 and 6. More pronounced lattice spacing effects are evident in the isentropic light-light and 
strange-strange quark number susceptibilities, which we attribute to the fact that these quantities 
have contributions only from the nonzeroth order Taylor expansion coefficients which are more 
sensitive to the cutoff. And finally, we did not find any peaks along the isentropic trajectories, 
which suggests that current experiments operate away from a possible critical point. 

A full-flavor quark-gluon plasma EOS is undoubtedly important for cosmological studies. Ac- 
cordingly, we determined the effects of the charm quark (at zero and nonzero chemical potential) 
and the bottom quark (at zero chemical potential only) on the EOS. Both heavy quarks were rep- 
resented in the heavy-quark-quenched approximation by asqtad valence quarks. We expect that 
the quenching error for such heavy quarks is small, especially for the b quark, but only a direct 
comparison with a calculation with dynamical c and b quarks can confirm that. We found that 
the contribution of the charm quark at zero chemical potential reaches about 20% in the energy 
density at temperatures of about 400 MeV and cannot be ignored in a high-precision cosmologi- 
cal calculation of the properties of the early Universe. The bottom quark contribution is within a 
standard deviation at that temperature. Our results for the charm and bottom effects on the EOS, 
however, may be affected by the heavy-quark discretization error we find in the free asqtad action 
calculation. This implies that they are possibly lower than their respective continuum values. 

At nonzero chemical potential, both charm and bottom quarks present a problem (the bottom 
quark much more so), since we found heavy-quark discretization effects in the Taylor expansion 
coefficients (especially large for the bottom quark), which could be overcome by tuning the coef- 
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ficient of the Naik term and/or using N t > 6. However, the charm and bottom quark contributions 
to the EOS due exclusively to the nonzero chemical potential are very small over the parameter 
range accessible to heavy-ion collisons, and at our level of precision they are entirely within the 
present statistical errors of the EOS at zero chemical potential. 
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Appendix A: General framework for adding the charm quark to the EOS at nonzero chemical 
potential 



With the addition of the charm quark to the u, d and s quarks in the sea, the partition function 
becomes: 

Z = j ^ e T lndet ^/ e T lndet ^v^ lndet ^~ 5 ^ (Al) 

with Mf being the quark matrix for flavor /. Thus, the pressure can be now expanded in the 
following manner: 

= 2^ c nmk{ 1 ) I j ; . j j ■ j 
1 n,m,k=0 \ / \ / \ 

where fif is the quark chemical potential for flavor / and the coefficients are 



n / - \ m / - \ k 



T 4 ~ 52 c nmk{T) ( j ( j ( j , 



(A2) 



Cnmk(T) 



1 1 1 N~ 



d n+m+k lnz 



(A3) 



with /jf is the quark chemical potential in lattice units. The coefficients above are nonzero only if 
n + m + k is even. A similar expansion applies to the interaction measure. These coefficients are 
for the asqtad quark action: 

' d{m f a) d n+m+k (2MJ l ) 



b, 



nmk 



1 N[ y 

s f=l,s,c 



~4 



dlna 



21 
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duo 



d\na 



tr 



d n+m+k( M -ldMl} 



Ml 



jc=0 d(piN t ) n d(n s N t ) m d(n c N t ) k 



1 N? 



d n+m+k (g) 



n\m\k\ Nf d(wN t )"d(v s N t ) m d(VcN t ) 



W,s,c=0 



In the above, £ = —dS g /d\na, with 5 g being the gluon part of the action. 



(A4) 



Appendix B: Calculating the pressure coefficients in the Taylor expansion 

These are most easily calculated using the following (similar to the derivation in the Appendix 
in Ref. [J): 



d/Ji 
dlnZ 

dlnZ 



A W0 = {U), 

-#oio = (Hi), 
*ooi = (fii). 



It can be shown that 



nmk 



dm 



■%n,m,k+l ■H-OOl-^-mnki 



where 

nmk - 

In all of the above, is defined as: 



e -Lo e -Ho e -Q 



Qk 



d n e L ° d m e H ° d k e Qo 

Brf Bjj" d/J k 

3 fc lndetM r 



(Bl) 
(B2) 
(B3) 



(B4) 
(B5) 
(B6) 

(B7) 
(B8) 



4 3/i* 

The operators and have a similar form for the light and the strange quark respectively. All 
coefficients c nm k which have at least one of the indices equal to zero have the same form as in the 
Appendix of Ref. [2], with appropriate substitutions of L n or H n with Q^. The "new" coefficients 
that appear to 0(6) are cn\ t \), C(3,2,i)> c 222 and C(4,i,i), where the notation (m,n,k) means all 
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distinct permutations of the indices. Explicitly we have 

C2ll = 2!1!1! jv^V^ 211 ~ --^011^200), (B9) 

° 222 = 2!2!2! A^V^^ 222 + 16 -^ 110 ^ 101 ^ 011 + Aa 101 #020 +4^oi i#200 + 4# n #002 + 2^002^200^020 

-^002^220 - -#200-#022 - -#202#020 - 4-#ioi-#121 ~ 4^110-^112 ~ 4^011-^211 ), (BIO) 

11 2 

C321 = 3!2!1 ,^ 3 ^ 3 (-#321 + 12^u ^ioi — 6^211-^110 — 2^011-^310 ~ 3^220-^101 ~ -#020-#301 

-3^200^121 +6^020-^200^ 101 + 12-#011-#110-#20o), (Bll) 
C411 = 4^777/^3 ^ 4U + 24 -^ HO-^101 -#200 - 4^301 ^110 -4^ 101 ^3 10 --#01 1-^400 

-6^200^211+ 6^01 i# 2 2 oo)- (B12) 

Permuting the indices above gives us the rest of the coefficients. Calculating the J? nm £ is straight- 
forward from Eq. (IB7I) . 



Appendix C: Calculating the interaction measure coefficients in the Taylor expansion 
1. First type of derivative 

This section gives a method to calculate the derivative: 



d n+m+k (M f l ) 



d{niN t ) n d(n s N t ) m d(n c N t ) k 



(CI) 



for f = l,s, c, when all of the indices n, m and k are nonzero. See Ref. [|2d for results when at least 
one is zero. It is convenient to define the observables: 



<B nmk = z-^- H ^-^^—L L———- \ (C2) 
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The above means: 



Sooo = (trM- 1 ), (C5) 
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<oo = (trM- 1 ). (C7) 
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Let / = / then we have the following rule: 
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Using the above, we calculate: 
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Replacing (B with <b' or <b" in the above, we get the expressions for the derivatives of (trM~ ) or 
(trA/^7 1 ). The explicit forms of f B n mk are easy to deduce from Eq. (IC2I ). To get the or 25^ 
we need to interchange appropriately the three observables: 
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along with L n , H n and <2« in the explicit forms of ( B n mk- 

2. Second type of derivative 

We also need to calculate the derivatives: 
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where again / = l,s,c. Similarly to the previous subsection, we define the observables 
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From the above, 
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Let / = /, then it is easy to see that 
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dCnmk 
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Similar expressions apply in the case of c' nmk and c" mk . Then the derivatives 
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have the form of the derivatives of (tr(M^ c ) S j in the previous section with the substitutions 
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3. Third type of derivative 

The last type of derivative that we need is the gauge derivative 

d n+m+k (g) 



d{wN t )»d{Li s N t ) m d{Li c N t ) k 

In this case, let 



(C32) 



W,s,c=0 



G „^ {ge -^ e -e>.——— } , (C33) 



and similarly as before 
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with 



Gooo = {§)• 

This means that the necessary derivatives ..A, ., u 
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have the same form as the 



derivatives 
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with 'Bnmk — > G nm k. The G nm ^ observables have very sim- 
ilar form to the A nm k observables, but with an additional multiplication by Q inside the ensemble 
average brackets of each term in them. 
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